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Abstract 

We calculate QCD correction to the neutralino annihilation cross section into quark 
anti-quark final state and discuss its implications to the calculation of neutralino relic 
density. We see that the QCD correction enhances the pair-annihilation cross sec- 
tion by O(10 %) when final-state quarks are non-relativistic. Consequently, when the 
lightest neutralinos dominantly annihilate into a ti pair, the relic density of the light- 
est neutralino is significantly affected by the QCD correction, in particular when the 
lightest-neutralino mass is close to the top-quark mass. 



1 Introduction 



Low energy supersymmetry has been regarded as an attractive candidate of new physics since 
it can solve various problems in the standard model of particle physics. One of the strong 
motivations of supersymmetry is that it provides a plausible candidate of cold dark matter. 
Recent precise cosmological observations have provided strong evidence of the existence of 
dark matter. In particular, Wilkinson Microwave Anisotropy Probe (WMAP) experiment 
has determined the present density parameter of the dark matter, which is given by pQ 



where hioo is Hubble constant in units of 100 km/sec/Mpc. From the particle-physics point 
of view, however, this raises a serious question since there is no viable candidate of dark 
matter in the particle content of the standard model. In supersymmetric models, the lightest 
superparticle (LSP) is stable due to -R-parity conservation and hence it can be a candidate of 
cold dark matter. It has been widely recognized that, if the lightest neutralino Xi is the LSP, 
its thermal relic density can become consistent with the WMAP value given in Eq. (jl.lj) in 
some parameter region. For a better understanding of the evolution of the universe, it is 
important to confirm or exclude the idea of neutralino (or LSP) dark matter. 

The first important step of such test is to find low-energy supersymmetry in future collider 
experiments. Once superparticles are found, in addition, it is also important to quantitatively 
test the idea of neutralino dark matter. Quantitative test will be performed by calculating the 
relic density of neutralino, using mass and coupling parameters of superparticles measured 
in collider experiments j2J. Currently, cosmological observations have determined the dark 
matter density with an accuracy of O(10 %), hence the density parameter of the neutralino 
Q x o is desired to be calculated theoretically at the same (or better) accuracy. 

As is well known, the relic density of neutralino is sensitive to the neutralino pair anni- 
hilation cross section, which strongly depends on mass and mixing parameters of the super- 
particles.* 1 Thus, for a precise calculation of the relic density, detailed calculation of the 
pair annihilation cross section is necessary. Indeed, in the case where the lightest neutralino 
is the LSP, the pair annihilation cross section and relic density of the lightest neutralino 
have been calculated in detail, in particular in connection with the WMAP observation [3]. 

So far, most of the previous calculations of pair annihilation cross section have been based 
on tree- level calculation. (See, however, 0j.) One may naively estimate the size of one- loop 
corrections to be 0(A 2 /167r 2 ) with A being relevant coupling constant, hence it is 0(1 %) or 
smaller as far as A ^ 1 . In such a case, radiative corrections to the pair annihilation cross 
section may be unimportant for the purpose of calculating Q x o with an accuracy of O(10 %). 
If we consider QCD correction, however, this is not always the case. If we consider processes 
with qq final state (with q and q being quark and anti-quark, respectively), QCD correction 
becomes O(10 %) when final-state quarks are non-relativistic. In the so-called "focus-point" 

^In some case, co-annihilation process also plays an important role. In this paper, however, we assume 
that the mass of the lightest neutralino is quite far from the second-to-the-lightest superparticle, and that 
the coannihilation process is ineffective. 
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model jS], for example, the annihilation process into a top quark pair, X1X1 ~^ ^, becomes 
one of the dominant pair annihilation processes [SJ U\. Then, QCD correction becomes 
important when the lightest neutralino mass is close to the top-quark mass due to threshold 
enhancement of the QCD correction. (The threshold enhancement is well known e.g. from 
the analysis of the cross section for e + e~ — > tt close to ti threshold [E]). 

The main purpose of this paper is to calculate the next-to-leading order (NLO) QCD 
correction to the neutralino annihilation cross section into qq final state and to discuss its 
implications to the calculation of the relic density of the lightest neutralino. (In fact, we also 
have to consider real gluon emission to take care of infrared divergences. Thus, we consider 
processes of the type X1X1 ~^ (wCsOn with g being the gluon.) As we will see, QCD correction 
may become O(10 %). Since the relic density of the lightest neutralino is approximately 
inversely proportional to pair annihilation cross section, this has certain significance for a 
precise calculation of the dark matter density. 

In addition, for scenarios where one of very weakly interacting superparticles (like grav- 
itino [9], axino ^U], or right-handed sneutrino ^T], which are recently called superWIMPs) 
is the LSP, study of the QCD correction has some significance. In such scenario, the light- 
est neutralino becomes unstable and decays into a superWIMP (and other standard-model 
particles), assuming that it is the lightest superparticle in the minimal-supersymmetric- 
standard-model (MSSM) sector. Then, cosmologically, the lightest neutralino becomes an 
important source of the superWIMP, which is now a candidate of cold dark matter. The 
important point is that the lifetime of the lightest neutralino is very long; in the thermal 
history of the universe, the decay occurs after the lightest neutralino freezes out from the 
thermal bath. In this case, abundance of the superWIMP produced by the decay of the 
lightest neutralino is given by the freeze-out density of the lightest neutralino. Thus, even 
if one of the superWIMPs is the LSP, precise determination of the freeze-out density of the 
lightest neutralino is also important.* 2 

This paper is organized as follows. In the next section, we summarize the framework 
of our study. Then, in Section El we present formulas for QCD correction to the process 
XiXi <?<?(<?)• 111 Section HI we calculate the relic density of the lightest neutralino and 
discuss implications of the QCD correction to the calculation of the dark matter density. 
Sectional is devoted to conclusions and discussion. 

2 Framework 

In our study, we work in the framework of the MSSM. In addition, we assume that the 
lightest neutralino xl becomes the LSP and hence is stable. In general, the mass and mixing 
parameters of superparticles depend on various MSSM parameters. To make our points 
clearer, we adopt several simplifications. 

In calculating fi x o, we need to fix properties of the lightest neutralino Xii which are 

# 2 One of the authors resisted glorifying supersymmetric models as scenarios for new physics, but the other 
authors ganged up and forced through the present style of Introduction. 
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determined from the neutralino mass matrix M. x o. Denoting the U(l)y and SU(2)l gaugino 
mass parameters as m G \ and m G 2 and supersymmetric Higgs mass as fin, M- x q * s gi ven by 



( —m G i — gifcos/3 g\vsva.(3 \ 

—Tn G 2 Q2V cos (3 —g2V sin /3 

—giVcos/3 g 2 vcos(3 /i# 

\ g±v sin /3 — sin (3 /ih / 



(2-1) 



where #i and c/ 2 are gauge coupling constants for for U(l)y and SU(2) L gauge groups, 
respectively, v ~ 174 GeV is the (total) vacuum expectation value of the Higgs boson, and 
tan (3 is the ratio of the vacuum expectation values of the up- and down-type Higgs bosons. 
This mass matrix is diagonalized by a unitary matrix U x o as 

U x oAi x oU x o = diag(m x o, m x o, m x o, m x o), (2.2) 

and the lightest neutralino is given by 

X? = [U x oY u B + [U xQ ]* 21 W° + [U x oY 31 H° d + [U^Hl (2.3) 

where B, W°, H® : and H® are Bino, (neutral) Wino, up-type Higgsino, and down-type 
Higgsino, respectively. 

In our study, we adopt the grand-unified-theory (GUT) relation among gaugino masses. 
In this case we obtain the relation between U(1)y and SU(2)l gaugino masses, m G i and 
m G2 , as 



m G 2 _ 3 m G i 
92 5 g\ 



(2.4) 



Then, since m G 2 > m G i, the lightest neutralino is approximately given by a linear combina- 
tion of Bino and Higgsinos. We will see that the QCD correction to the process XiXi ~^ ^G?) 
becomes most significant when the lightest neutralino is mostly Bino-like but is with a siz- 
able Higgsino contamination. Notice that the formulas for the QCD correction which will 
be given in the next section are independent of the relation among the gauginos. 

With the unitary matrix U x o, the Z~x\~x\ vertex, which is relevant for the process XiXi ~^ 
qq(g), is given in the form 



1 

2' 



^x°x° = n A xXil^5XiZf„ (2.5) 



where 



A x = \g z ([U x „r 31 [U x o} 31 - [U^ [U x o} 41 ) , (2.6) 



with g z = yj g\ + g\- When m G1 <C fi H , [U x o] u ~ 1 while [£/ x °]ii (i = 2,3,4) become 
close to 0. In this case, the lightest neutralino is almost Bino and its interactions with 
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gauge and Higgs bosons are suppressed. Consequently, dominant annihilation processes are 
into lepton pairs (as far as the sleptons are lighter than the squarks) and hence are p-wave 
processes. When \xn is relatively small, on the contrary, Higgsino component in the lightest 
neutralino becomes enhanced and [£/ x °]3i and [t^ x o] 4 i become sizable. In this case, the lightest 
neutralino pair dominantly annihilates into ti, W + W~, and ZZ final states since they are 
s-wave processes. 

We also present the interaction of the neutralino with the lightest (standard-model-like) 
Higgs boson h. The interaction term is given in the form 

£h X o x o = \y x x\x% (2.7) 

In our following study, we consider the so-called decoupling limit where all the physical 
Higgs bosons except the standard-model-like Higgs h are very heavy. In this case, we can 
neglect effects of heavier Higgses on the annihilation processes of the lightest neutralino. In 
addition, in the decoupling limit, Y x is well approximated by* 3 

Yx = ~ {92 [U x o} 21 - 9l [U x o] n ) ([U x o} 31 cos/5 - [C/ x o] 41 sin/?) . (2.8) 

(Here and hereafter, we work in the bases where the matrix U x o is real.) 

In the calculation of f^o, sfermions may also contribute. In the following, we pay par- 
ticular attention to the case where the annihilation process of Xi is dominated by s-wave 
processes. In such a case, sfermion-exchange diagrams are mostly subdominant since they 
induce p-wave processes. One exception may be t-channel exchange of stops. We assume 
that the stops are heavy and their contribution to the pair annihilation of Xi is negligible. 
In addition, we also assume that heavier Higgs bosons are heavy enough so that they do not 
significantly affect Q x o . Such a mass spectrum with heavy scalars is realized in, for example, 
focus-point models 



3 QCD Correction to 0" x o x o_> g ^) 

We calculate the NLO QCD correction to the cross section for neutralino pair annihilation 
into qq(g) final state. In this paper, we consider the case where the squarks (in particular, 
stops) are much heavier than the Z-boson. This is the case in a large class of supersymmetric 
models, partially due to renormalization group effects via gluino mass. In such models, the 
decays X1X1 ~* QQ{9) are induced dominantly via Z and /i-boson exchange diagrams. When 
the initial-state neutralinos are non-relativistic, Z-boson exchange diagram turns out to be 
most important. Z and /i-boson exchange contributions do not interfere, since the P (parity) 
properties of Z and h are different; hence we may discuss their contributions separately. 

# 3 In the general case, Y x is given by replacing j3 — > a + ^ir, where a is the mixing angle in the CP-even 
Higgs sector. (See, for example, [T^].1 
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3.1 0(a s ) correction 

We first consider 0(a s ) contributions from the one- loop and real-gluon emission processes. 
As will be discussed below, QCD corrections are enhanced when the final-state quarks are 
non-relativistic due to boundstate effects. In such a case, it is necessary to take account of 
(some part of) higher-order corrections in order to obtain reliable predictions. This will be 
treated in the next subsection. 

We start with the contribution of the Z-boson exchange diagrams. The relevant interac- 
tion of the Z-boson with quarks can be parameterized as 

C-Zqq = J2 + A <ll$)<l Z ' p- (3- 1 ) 

q 

As we will see, in the calculation of the relic density of Xv tt(g) final state is of particular 
interest. For the top quark, 

Vt= w£ A -* A )' A=A * sz - (3 ' 2) 

We should consider two contributions simultaneously; one is the virtual correction to Z- 
q-q vertex and the other is the real gluon emission process XiXi ~ > Z* ~^ Both of them, 
individually, are infrared divergent, but these divergences cancel when both contributions 
are summed. Let us denote the Feynman amplitude for the process XiXi ~ ¥ Z* ~ > <?<?(<?) as 

M x ^ z ^ {g) =A^D^Aft. (3.3) 

Here, Dffl = Q 2 ^ m '2 (g^iu — ) is the propagator of Z-boson (with being the total four- 
momentum of the system), while Aft) and Aft are bi-spinors (including coupling constants) 
consisting of the wave functions of the neutralino and quark, respectively. It is convenient 
to define the following quantity: 

tt(?) V Id®* Ia (<?) A (i? H + V ( d<& - lA (g) A (9 H (3 A) 



UAl 

spin spin 



where the first and second terms represent contributions from the processes XiXi Z* ~~ > 
qq (i.e., leading order + virtual correction) and XiXi ~^ Z* — >• qqg (i.e., real emission), 
respectively. In the above formula, the integrals J c?$ are carried out over the final-state phase 
space, and the spin sums of the final-state particles are taken. The first and second terms in 
Eq. ()3.4|) individually contain infrared divergences, while lift is finite. We decompose lift 
into the vector and scalar parts as 

n<g = n^(Q 2 ) (-g, v + 9&l\ + Tlft(Qi)^. (3.5) 
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(a) 



(b) 




For the purpose of the calculation of the relic density of neutralino, we define* 4 



2 V xl a xWi^ Z *^(9) 



4Q 2 



spin spin 



-'-IQi-l 



(3.6) 



where v x o is the velocity of the lightest neutralino in the center-of-mass frame, which is given 
by 



v x1 



# 4 It is often the case that the quantity S is expressed as v IC \a, where u re i is the relative velocity of the 
initial-state neutralinos and a is the neutralino pair-annihilation cross section. 



Am? o 



(3.7) 
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In addition, <J x o x o 



z*^qq( g ) is the total cross section for the process XiX° ~~ * Z* — > qq(g), 
hence, in the final formula, the spins of neutralinos are understood to be averaged. As we 

z*~*qq( g ) plays an important role in the calculation of f2 x o. With 

is expressed as 



will see, the function S x o x o 



Ilty and II S 



(<z) 



S x o x o_ 



S x o x o_ 



*qq{g) 



(Q 2 



1 

Q 1 



Q 2 



Am 2 



(?) 



(Q 2 



m 



2 ^2 ii v 



(Q 2 ) + 



2m 2 



VTLr 



14 



(3. 



Importantly, the first term in the parenthesis is proportional to Q 2 — 4m 2 ; as a result, it 



is suppressed when the initial state neutralinos are non-relativistic. Thus, in such a case, 



the contribution of ITy' ) becomes negligible. In addition, Ilg''' is proportional to the mass- 
squared of the final-state quark. Since the relic density is determined by the cross section 
at Q 2 ~ Am 2 , we pay particular attention to the process X1X1 ~ ¥ Z* ~ ¥ tt(g) an d neglect 
other processes concerning lighter quarks, in calculating the relic density of Xi- 

In this subsection, we present the formulas for Ily^ and II^ up to 0(a s ). Corresponding 
Feynman diagrams are shown in Fig. ^ At the tree level, we obtain* 5 

,2 



(«) 



U^(Q 2 



tree 



rt«\Q 2 



3tt(1 

2 



(3-«J)|V s 



\ 2 + 2v 2 \A g \ 2 



m 



71 



- V q \ A q\ 



where 



4m 2 



\ Q 2 



(3.9) 
(3.10) 

(3.11) 



All information on the QCD correction (at 0(a s )) is contained in Ily and ITg^. In the 
calculation of the virtual contribution (Fig.^b)), we adopt the on- mass-shell renormalization 
condition to determine the vertex counter term: we require that the 0(a s ) correction to the 
vector vertex (which is proportional to V q ) vanishes when Q 2 —>■ 0. Then, 0(a s ) correction 
to IT^ is found to be [T3| 



0(a s ) 



a s C 2 m 2 q 2 

2^ \ V 1\ 



24vr 2 (l 



Vi 



1(3 -v 2 )A(v q ) + 6v q (5 -3v 2 ) 



2 - 7 ..4^ 1 + v q 



-(33 + 22^-7^) In 



1-v, 



q J 



asC 2 m 2 q 2 



48vr 2 (l 

+ (21 + 59v 2 + Wv* - 3v*) In 



32v 2 A(v q ) + Qv q (-7 + 10v 2 + v 4 q ) 

l + Va 



(3.12) 



# 5 At 0(a s ), the correction from real gluon emission is also from tree diagrams. [Ily ] tree an d [Ilg ] tree do 
not include such contributions. 
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where C2 = t is quadratic Casimir operator for the fundamental representation of SU(3), 



A(v q ) 



'X + vi 



4Li, 



1 — V n 



l + Va 



+ 2Li 2 



-2 In v „ In 



l+v g 



3v q In — 



1 + V a J l + V 1—Va 



— - 4w g lnw 9 , 



and Li2(x) is the dilogarithm (Spence) function 

Li 2 (V) = - 



* ln(l-i) 
dt— '-. 



(3.13) 



(3.14) 



In addition, 0(a s ) correction to Ilg^ is given by 



0(a 3 ) 



a s C 2 m 2 q 

167T 2 



\A n 



l§A{v q ) + Qv q {7 - v 2 ) + (19 + 2v 2 q + 3v 4 q ) In — — * 

1 Vn 



(3.15) 



We have checked that the ratio [H^\Q 2 )}o(a B )/[^s > (Q 2 )]tree coincides with the correction 
factor for the decay rate of the pseudo-scalar Higgs boson given in [T4] . 

We can also calculate the /i-exchange contributions. Although the process ~~ * h* — > 
ti(g) is p-wave suppressed and hence is subdominant, we write down the formula for the 
QCD correction to this process to make this paper self-contained. Let us denote the hqq 
vertex as 



(3.16) 



and define 



S x o x o^ h ,^ qq{g) (Q 2 ) - 2v x oa x o x o^ h „^ qq{g) . 



(3.17) 



Adopting on-mass-shell renormalization condition for the h-q-q vertex, the coupling constant 
Y q is related to the (on-shell) mass of q as 



Yn 



V2v' 



Then, S x o x o_ h *_> qq - {g) is given by 



S x o x o^ h ,^ qq{g) 



m 2 (Q 2 - 4m V 



(O 2 ) = ~'" Xl v 3 \YY I 

W j 8nQ 2 (Q 2 -m 2 ) 2V(il q xl 



A, 



71 



where 



A, = U Vq ) + * (3 + 34^ - 13<) In + * (-1 + 7l fl 



(3.18) 



(3.19) 



(3.20) 



For the calculation of the relic density of neutralino, the s-wave contribution plays the 
most significant role. As we mentioned, for the s-wave process, ti(g) final state is most 
important while the effects of lighter quarks are unimportant (if we restrict ourselves to the 
Z-exchange diagrams). Thus, in order to see the size of the QCD correction to the s-wave 
part, we define 

P (r>2\ _ [ n s } (Q 2 )]o(q s ) /oo-n 

Rot ° JQ } " tnJ>(o>)u ' (3 ' 

We can see that the 0(a s ) correction to the s-wave part is enhanced when vt becomes smaller. 
Indeed, expanding Ro(a s ) around v t = 0, we obtain 



D C 2 a* 
R-o(a 3 ) ■ 



7T 



TT 2 



2vt 3 + OK) 



(3.22) 



This ratio is inversely proportional to Vt when vt < 1. As a result, the QCD correction 
largely enhances the neutralino annihilation cross section especially when the final-state top 
and anti-top quarks are non-relativistic. Thus, in the calculation of the relic density of the 
lightest neutralino, QCD correction may become important. 

In fact, it is well known that QCD corrections are enhanced as v t becomes smaller, due 
to boundstate effects. Naive perturbative expansion breaks down in the region v t ^5 ot s . In 
this region, systematic treatment of boundstate effects are known for heavy quarks, such as 
top quark. Up to NLO, this corresponds (formally) to the resummation of terms of the form 
(a s /v t ) k and v t (a s /v t ) k = a s (a ;s /f t ) fc ~ 1 for < k < oo (where k is an integer). In the next 
subsection, we compute the QCD correction reliably by resumming all these corrections. 



3.2 NLO annihilation cross section close to qq threshold 

We compute the QCD correction to the process XiXi — > <?<?(<?) close to qq threshold up to 
NLO (including boundstate effects). Our main concern is to calculate the relic density of the 
lightest neutralino; for this purpose, the s-wave part of the cross section is most important. 
Thus, in this subsection, we concentrate on the QCD correction to the s-wave part. 

Following the standard framework developed for the threshold cross sections, for example, 
for the process e + e~ — > ti jHj, the NLO annihilation cross section (when the final-state quarks 
are non-relativistic) is given as follows: 



NLO, NR 



1 + 



C 2 a s 



A 



(hard) 



7T 



(3.23) 



(hard) 

where A s represents the hard-vertex correction given by 



(hard) 



lim 



C 2 a, 



TT 



-1 



n^(Q 2 ) 



\0(a s ) 



J tree 



7T 

2% 



-3. 



(3.24) 
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Boundstate effects are included in F(Q 2 ), which is given in terms of the Green functions of 
the non-relativistic Schrodinger equation: 



Im[G(0; E)} 
Im[G (0;E)]' 



G(x; E) is defined by 



(E + iV q 



V 2 



+ VQCD (^Bohr) 



G(x; E) = 5\x). 



(3.25) 



(3.26) 



Here, Vq^d (r; /iBohr) denotes the static QCD potential up to NLO 15j E = — 2m q is 
the center-of-mass energy measured from the qq threshold, and Y q is the total decay width 
of q. On the other hand, Gq(x; E) is the non-relativistic Green function of a free qq pair, 
which is defined via Eq. (j3~27^l after setting V^ 0) (r) to 0. 

In the case of QED (i.e., when Vqcd(^) is replaced by the Coulomb potential — c^qed/?", 
with «qed being the fine-structure constant), F(Q 2 ) reduces to the well-known Sommerfeld 
factor: 



QED 



1-e" 



where 



vr «QED 



(3.27) 



(3.28) 



In this case, it is easy to see via Eqs. ()3.22|) #6 and (|3.27|) that the QED correction including 
the boundstate effects is twice as large as the (naive) one-loop correction at \[Cp = 2m q . 
One may regard this as a reference for the significance of boundstate effects, with respect to 
the naive one-loop correction, in the limit where the running of the strong coupling constant 
may be neglected (since it corresponds to the limit where the static QCD potential reduces to 
the Coulomb potential); qualitatively it is a reasonable approximation for final-state quarks 
as heavy as the top quark. 

With Eq. (|3.23|) . we also define the magnitude of the QCD correction relative to the LO 
contribution as 



R(Q 2 ) NLO, NR 



I NLO, NR 



[i4 9) (q 2 )u 



- 1. 



^3.29) 



3.3 QCD-corrected neutralino annihilation cross section 

In calculating the relic density of the lightest neutralino, the annihilation cross section for 
Q 2 ~ 4m 2 o is important. This is because the lightest neutralinos decouple from the ther- 
mal bath when the cosmic temperature is much lower than m x o; at such low temperature, 
neutralinos are non-relativistic. 

# 6 One should replace Cia s by q;qed- 
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To see the importance of the QCD correction in such a case, in Fig. El we plot Ro(a a ){Q 2 = 
Am\) and -Rnlo nr(Q 2 = 4m^ ) as functions of m o. It can be easily seen that the QCD 
correction significantly enhances the cross section. Since the tree-level cross section vanishes 
as v t — > while the QCD corrections stay constant, Ro{a a ){Q 2 = 4m^ ) and -Rnlo,nr(Q 2 = 
4m^ ) diverge as m x o — > m t . As a result, the size of the QCD correction becomes comparable 
to the tree-level contribution when the final-state top and anti-top quarks are non-relativistic. 
This fact has important implications to the calculation of Q x o when the lightest neutralino 
dominantly annihilates into a tt pair. 

Estimate of the annihilation cross section using the one-loop QCD correction Ro(a s ) is 
valid when a s <C Vt and a 8 <l, since both (naive) higher-order contributions a™ and higher- 
order boundstate contributions (a s /v t ) n are suppressed. On the other hand, estimate of the 
QCD correction by -Rnlo.nr is valid when v t , a s <C 1. Thus, in the overlap region when 
a s < tit < 1 hold, both Ro(a s ) an d -Rnlo,nr are reasonable approximations. According 
to Fig. El when the lightest neutralino mass is in the range 200 — 220 GeV, Ro(a a ){Q 2 = 
Am 2 ) is reasonably close to -Rnlo nr(Q 2 = 4m^, ) (with respect to the estimate of the 
NNLO correction). On the contrary, with smaller m x o, boundstate effects become important. 
Indeed, in such a case, Ro(a a )(Q 2 = 4m^ ) and -Rnlo nr(<5 2 = 4m^ () ) differ significantly, and 
the latter prediction is more reliable. 

Since the QCD correction enhances the cross section & x xl-*ti(g)> it ma y pl a Y an important 
role in the calculation of the relic density of Xi- in the next section, we will discuss how the 
relic density changes as we take account of the QCD correction. 



4 Relic Density of the Lightest Neutralino 

In this section, we calculate the relic density of the lightest neutralino. We include the QCD 
correction to the neutralino annihilation cross section given in the previous section and see 
its effects on the calculation of the relic density of the neutralino. 

QCD corrections for the process XiXi ^ QQio) ma ¥ be important for the precise cal- 
culation of f2 x o. When the MSSM parameters will be precisely measured in future collider 
experiments after the discoveries of superparticles, cross sections of the pair annihilation 
processes of \i w iU be calculated in detail. Using the formulas given in the previous section, 
it is straightforward to include the QCD correction to the processes XiXi ~" * QQ(q)- Here, we 
show the importance of the QCD correction. To make our point clearer, here we take into 
account only the most important annihilation processes and neglect sub-dominant ones. 

In order to calculate the relic density, we should follow the evolution of the number density 
of the lightest neutralino n x o, which is governed by the Boltzmann equation. Denoting the 
energy distribution function of the lightest neutralino as / x o,# 7 the Boltzmann equation is 

d 3 p o 

# 7 n x o and f x o are related to each other by n x o = J ( 27r )3 /%g (Px° ) ' 
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Figure 2: QCD correction to the cross section ^ x \x\-^m{g) wnen the initial \i are at res t (i- e -> 
Q 2 — > 4m^ ) as a function of the LSP mass m x o. Three lines correspond to Ro{a s ) (solid), 
-Rnlo.nr with /iBohr = 20 GeV (dashed), and -Rnlo,nr with /iBohr = 40 GeV (dotted). 
The input values are taken as a s (m t ) = 0.1080 (corresponding to a s {mz) = 0.1187) and 
m t = 172.5 GeV. /^Bohr represents the renormalization scale in VQ^\r), where a sensible 
choice of /^Bohr is at the Bohr scale of ti resonances. We may regard the difference between 
the predictions for /^Bohr = 20 GeV and 40 GeV as an estimate of contributions beyond NLO. 



given 



by jH] 



dn^o 



processes 



X-l Xi 



+ (production term), 



(4.1) 



where H is the expansion rate of the universe. In the above formula, / o?$f is the integration 
over the phase space of final-state particles and we sum over all the relevant pair annihilation 
processes. In addition, the phase space integrations / c/ILo and / dU' ' are for the initial-state 
neutralinos (whose momenta are p x o and p' x o, respectively). 

Even if the superparticles may not be in chemical equilibrium, they are expected to be 
in kinematical equilibrium since superparticles effectively interact with ordinary standard- 
model particles in thermal bath. In addition, since the freeze out of the lightest neutralino 
occurs when T <C m v o, we use the Boltzmann distribution function for neutralinos; at cosmic 

Xi ' 
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temperature T, 



/ x j( Px o) = 2e-H--^, (4.2) 



where the factor of 2 represents the spin multiplicity of neutralino, p x o is the energy of x%, 
and n x o is the chemical potential. 

In order to follow the evolution of n v o for T <C m v o, it is convenient to define 

x(T) = — . (4.3) 
m v o 

Once the temperature becomes lower than m x o, x can be used as an expansion parameter; 
we expand various quantities using x. For example, the number density of the lightest 
neutralino is given by 



n^o = 2e 



/ 3/2 
-(m x O-M x o)/T /m x? M 



1 - —x + 0{x 2 



(4.4) 



where this formula is obtained by expanding the energy as 

\Px°f \Pxl\ 4 

p y o o = m Y o + — — l - 1 H , 

^xi.u xi 2m Y o 8m 3 

with p x o being the three-momentum of x?- Substituting Eq. ()4.2j) into Eq. ([4.1)1 . and applying 
the detailed balance theorem to determine the production term, the Boltzmann equation 
reduces to 

+ 3ifn x? = " (l - ^) 5 tot (g 2 = 4m x? + 6xm x? ) (n^ - n^ )2 ) + 0(x 2 ), (4.5) 

where is the chemical-equilibrium value of the number density of Xi: which is obtained 
from Eq. (|4.4jl by taking /i x o — > 0. In addition, «S to t is given by 

Stot = 2v x a x o x o^ aU , (4.6) 



with a' x o x o_ >all being the total pair annihilation cross section of the lightest neutralino.* 8 
Taking account of the correction of the order of x, Q 2 is shifted from 4m 2 to 4m 2 + 6xm 2 . 
This is due to the fact that, in the thermal bath, neutralinos have non- vanishing kinetic 
energy whose average is |T. 



# 8 Here, we assume that the effects of coannihilation is negligible. For a precise calculation of the relic 
density, effects of the coannihilation should be taken into account ^7|. Once the relevant MSSM parameters 
are measured at the future collider experiments, however, effects of the coannihilation can be theoretically 
controlled. Thus, in our analysis, we do not include effects of the coannihilation, which are beyond our scope. 
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From Eq. (|4.5|) . several important features of n x o are derived. When the cosmic tem- 
perature is high enough, production and annihilation rates of the lightest neutralino, which 
are ~ St a tn x o, are relatively high. At this epoch, scattering rate becomes much larger than 
the cosmic expansion rate H and, consequently, the number density of Xi keeps up with 



the equilibrium value n^. As the temperature decreases, however, the situation changes. 

,(eq) 



Importantly, once the cosmic temperature becomes lower than m x o, n^o is Boltzmann sup- 
pressed. When n y o (which is close to n o at this epoch) becomes so suppressed that S tot n x o 
becomes smaller than H, Xi freezes out from the thermal bath. After the freeze out, the 
number of Xi m comoving volume is conserved. 

The freeze-out of the lightest neutralino occurs when x _1 ~ 20 and hence the lightest 
neutralinos are non-relativistic in thermal bath at the time of freeze out. Thus, the relic 
density is primarily determined by the behavior of the total cross section at Q 2 — > 4m^o, as 
can be understood from Eq. (|4.5J) . When A x given in Eq. ()2.6|) is sizable, pair annihilation 
of the lightest neutralino is dominated by s-wave processes as Q 2 — > 4m 2 . Thus, in such 

Xi 

a case, il x o is mostly determined by the s-wave part of the total cross section. When all 
the sfermions and heavier Higgses are relatively heavy, the relevant s-wave processes are 
XiXi ~~ * tt(g)i W + W~ , and ZZ. In our analysis, we assume that this is the case and 
approximate 

5 tot ^ S xlx\-*Z*-+qq{ g ) + + ^"x^ZZ- (4-7) 

Here, the second and third terms on the right-hand side represent the contributions from 
W + W~ and ZZ final states, respectively: 

^x° 1 Xi^ w+w ~' zz = ^x^xlxl^w+w-^zi (4-8) 

with cr x » x o^w+w-,zz being total cross section for the corresponding processes. We use the 
tree-level formulas for a x o x o_> w+w - and a x o x o_, zz . In addition, we evaluate these quantities 
at the threshold since their Q 2 dependences are weak when x C 1. (For XiXi ~ > ^(sO? 0(x) 
terms are considered unless otherwise mentioned.) Other processes, like XiXi ~^ Zh and hh, 
are p-wave processes and are neglected in our study. If we include those processes, Q x o may 
decrease by a few % or so [T7j; once the MSSM parameters will be experimentally determined, 
it will be straightforward to take account of effects of the p-wave processes. Furthermore, 
Sx°Xi^qq(g) ^ s calculated neglecting the /i-exchange diagram since its contribution is p-wave 
suppressed. 

We have numerically solved the Boltzmann equation ()4.5|) to determine the present value 
of n x o. We impose the initial condition n x o = n^f* at the cosmic time where «StotW^o ^> H 
and follow the evolution of n x o. Then, we determine the number density of the lightest 
neutralino at S to tn^ <S H. In our following argument, it is convenient to define the yield 

Xi 

variable 

Y x o = (4.9) 
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with s being the entropy density. Important property of this quantity is that Y x o stays 
constant when S to t n ^ ^ H: indeed, we have numerically checked that the yield variable 

Xi 

becomes a constant of time when the scattering rate is much smaller than the expansion 
rate. Using the yield variable after the freeze-out of Xi (which we denote Y^™^), the 
density parameter is given by 



T/(now) 



Pcrit 1 1 



(4.10) 



t (now) 

Here, p cr it is the critical density and s( now ) is the present entropy density; numerically, 

Pcrit - 3.6 x 10- 9 /4o GeV. (4.11) 



-.(now) 



In Fig. El we plot the density parameter of the lightest neutralino as a function of 



m v o. 

In the figure, we compare the results with and without the QCD correction to the process 
XiXi ~^ ^"(flOi we evaluate the cross section & x ° x -+ti(g) with NLO QCD correction with 
boundstate effects (solid), 0(a s ) correction (dashed), and tree-level formula (dotted). In the 
calculation with the boundstate effects, we restrict ourselves to the case with m x o ^ 210 GeV. 
This is because, when m x o is too large, non-relativistic approximation used in the calculation 
becomes unreliable. Fortunately, however, as the lightest neutralino mass becomes larger 
than ~ 200 GeV, result with the boundstate effects is in a very good agreement with that 
with 0(a s ) contributions. This fact shows that QCD correction to the process XiXi — > tt(g) 
is well approximated by the one-loop formula when the lightest neutralino mass is larger 
than ~ 200 GeV. On the contrary, when m x o is smaller, difference between the one-loop 
result and the resummed one is quite large; in this case the latter prediction (solid) is more 
reliable. 

We can see that the QCD correction significantly reduces the density parameter of the 
lightest neutralino when m x o is close to the top-quark mass. Since the ti final state is more 
important than the gauge-boson final states when the lightest neutralino mass is relatively 
close to the top-quark mass, this fact can be easily understood from Fig. |2J As m x o becomes 
larger, on the contrary, <r x o x o^w+w-,zz becomes larger than ff x ° x °^tf(g) an d the gauge-boson 
final state becomes more important than the ti final state. In such a case, even if the QCD 
correction to the process XiXi ~~ * tt(g) is sizable, Q x o may not be affected so much. Indeed, 
we can see such a behavior in the figure. As the lightest neutralino mass becomes large 
enough, tree-level result becomes almost the same as the one with QCD correction although 
the QCD correction to a x o x o^ t i^ is still more than 10 %. 

In Fig. Efa), we kept 0(x) terms of S x o x o^ t ^ g y For comparison, we also calculated the 
density parameter with neglecting x-dependences (by taking x — > 0); the results are shown 
in Fig. ISJb). We can see that the resultant density parameters differ in two figures when 
the lightest neutralino mass becomes close to the top-quark mass. We have checked that the 
effects of 0(x) terms become larger as m x o becomes closer to m t . 

In Figs. |U and we also plot the density parameter for the cases with fi H = 260 GeV 
and ix h = 300 GeV, respectively. As we increase hh, the WMAP value is realized with 
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Figure 3: Density parameter of the lightest neutralino as a function of rn x o. The cross 
section a x o x o_^ t ^ is evaluated with NLO QCD correction with boundstate effects (solid), 
0(a s ) correction (dashed), and tree-level formula (dotted). Here we take a s {m z ) = 0.1189, 
m t = 172.5 GeV, and fin — 280 GeV. The dark shaded band shows the dark matter density 
measured by the WMAP (see Eq. (jl.ljl ). For (a) and (b), collision term for ti(g) final state 
is evaluated up to 0(x) and O(x ), respectively. 
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larger value of m x o due to the change of the mixing matrix of the neutralinos. These figures 
show that, with these values of the relic density of the lightest neutralino is significantly 
affected when m x o is close to m t . As hh increases, however, QCD correction becomes less 
important in the parameter region where Q x o agrees with the WMAP value. 

5 Conclusions and Discussion 

In this paper, we have calculated QCD correction to the process X1X1 ~^ Q<l(9)i an d discussed 
its implications to the determination of fl x o. We have seen that the QCD correction enhances 
the cross section o~ x o x o_ O(10 %) when the final-state quarks are non-relativistic. 
Importantly, the lightest neutralino dominantly annihilates into a tt pair (as well as into 
W + W~ and ZZ) when the lightest neutralino has sizable Higgsino components. In such 
a case, QCD correction may become large and significantly enhance the total annihilation 
cross section of the lightest neutralino. Since the QCD correction is enhanced close to the 
threshold by factors (a s /v t ) k , the QCD correction becomes important in particular when 
m x o is close to the top-quark mass. Since the relic density of the lightest neutralino is 
(approximately) inversely proportional to the total pair annihilation cross section, the relic 
density may be reduced by 0(10 %) by the QCD correction. Thus, for a precise calculation 
of the relic density of the LSP, QCD correction should be taken into account when Xi 
dominantly annihilates into the tt final state. 

Since the dark matter density has been determined with an accuracy of 0(10 %) by 
WMAP and other experiments, this fact has important implications for future studies of 
supersymmetric dark matter (s). Once properties of superparticles are determined in future 
collider experiments, it will become possible to theoretically calculate the relic density of 
the lightest neutralino (if it is the LSP). Comparison of the result of such a calculation 
with the observed density of the dark matter will be a crucial test of the scenario of the 
lightest-neutralino dark matter. In addition, even if the lightest neutralino is not the LSP, 
precise determination of the freeze-out density of the lightest neutralino may be important 
for understanding the origin of the dark matter. In particular, possibilities of very weakly 
interacting LSP (like gravitino-, axino-, or right-handed-sneutrino-LSP), which are called 
superWIMPs, have been intensively considered recently. If one of the superWIMPs is the 
LSP, the lightest neutralino becomes unstable. In this case, however, the lifetime of the 
lightest neutralino is usually very long and hence, if Xi is the lightest superparticle in the 
MSSM sector, thermally produced lightest neutralinos decay after they freeze out from the 
thermal bath. If this is the case, thermally produced lightest neutralino becomes one of the 
significant sources of the superWIMP cold dark matter. Thus, even in the scenario with 
superWIMPs, it is important to calculate the relic density of the lightest neutralino. 

Acknowledgements: The authors would like to thank T. Asaka, K. Hikasa and M. Tanabashi 
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